//////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////// 
// Create a base axis on boundary of surface mesh

#include "baseOnBoundaryMesh.h"
#include "geom/point3d.h"

//////////////////////////////////////////////////////////////////////////
baseOnBoundaryMesh::baseOnBoundaryMesh(void)
{
}


//////////////////////////////////////////////////////////////////////////
baseOnBoundaryMesh::~baseOnBoundaryMesh(void)
{
}

//////////////////////////////////////////////////////////////////////////
bool baseOnBoundaryMesh::process( cTopology *top )
{
  int i, idi, idj, edge;

  // check for number of loops
  int nloops = top->GetNumLoops();
  if (nloops != 1)
    return false;

  // get loops 
  int npts, *id_pts;
  top->GetLoop(0, &npts, &id_pts);

  // get all points
  double coord[3];
  Point3D *pts = new Point3D[npts];
  for (i = 0; i < npts; ++i)
  {
    top->GetCoordNode(id_pts[i], coord);
    pts[i].setXYZ(coord[0], coord[1], coord[2]);
  }

  // resize all vectors
  m_edges.resize(npts);
  m_edges.resize(npts);
  m_idEdges.resize(npts);
  m_idNodes.resize(npts);

  Point3D vec_edge_elem, vecij, edge_center, tmp;

  // local axis for edges
  for (i = 0; i < npts; ++i)
  {
    idi = id_pts[i];
    idj = id_pts[(i+1)%npts];
    edge = top->GetEdge(idi, idj);

    // set id edge
    m_idEdges[i] = edge;

    // vector in the edge from i to j
    vecij = pts[(i+1)%npts] - pts[i];
    edge_center = (pts[(i+1)%npts] + pts[i]) * 0.5;

    // adjacent element to edge
    int adj_faces = top->AdjElemToEdge (edge, 0);

    // center of face
    top->GetElemCenter (adj_faces, coord);
    Point3D elem_center (coord[0], coord[1], coord[2]);

    // vector form edge to center of face
    vec_edge_elem = elem_center - edge_center;

    // set local axis w, direction of edge
    m_edges[i].w = vecij;
    m_edges[i].w.normalize();

    // set local axis v, normal to face
    m_edges[i].v = vec_edge_elem.crossProd(vecij);
    m_edges[i].v.normalize();

    // set local axis u, tangent to face in out direction
    m_edges[i].u = m_edges[i].v.crossProd(m_edges[i].w);

    // center of transformation
    m_edges[i].center = edge_center;
  }

  // local axis for nodes - average
  for (i = 0; i < npts; ++i)
  {
    m_idNodes[i] = id_pts[i];

    m_nodes[i].center = pts[i];

    m_nodes[i].u = (m_edges[i].u + m_edges[(i+1)%npts].u) * 0.5;
    m_nodes[i].v = (m_edges[i].v + m_edges[(i+1)%npts].v) * 0.5;
    m_nodes[i].w = (m_edges[i].w + m_edges[(i+1)%npts].w) * 0.5;
  }

  delete []pts;

  return true;
}
